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ABSTRACT 

We combine the six high-resolution Aquarius dark matter simulations with a semi-analytic 
galaxy formation model to investigate the properties of the satellites of Milky Way-like 
galaxies. We hnd good correspondence with the observed luminosity function, luminosity- 
metallicity relation and radial distribution of the Milky Way satellites. The star formation 
histories of the dwarf galaxies in our model vary widely, in accordance with what is seen ob- 
servationally. Some systems are dominated by old populations, whereas others are dominated 
by intermediate populations; star formation histories can either be continuous or more bursty. 
Ram-pressure stripping of hot gas from the satellites leaves a clear imprint of the environment 
on the characteristics of a dwarf galaxy. We hnd that the fraction of satellites dominated by old 
populations of stars matches observations well. However, the internal metallicity distributions 
of the model satellites appear to be narrower than observed. This may indicate limitations in 
our treatment of chemical enrichment, which is based on the instantaneous recycling approxi- 
mation. We hnd a strong correlation between the number of satellites and the dark matter mass 
of the host halo. Our model works best if the dark matter halo of the Milky Way has a mass of 
~ 8 X lO^^M©, in agreement with the lower estimates from observations, but about a factor of 
two lower than estimates based on the Local Group timing argument or abundance matching 
techniques. The galaxy that resembles the Milky Way the most also has the best matching 
satellite luminosity function, although it does not contain an object as bright as the SMC or 
LMC. Compared to other semi-analytic models and abundance matching relations we find 
that central galaxies reside in less massive haloes, but the halo mass-stellar mass relation in 
our model is consistent both with hydrodynamical simulations and with recent observations. 

Key words: galaxies: evolution - galaxies: formation - galaxies: dwarf - galaxies: abun- 
dances - galaxies: stellar content - Galaxy: halo 



1 INTRODUCTION 

There is much to be learned about galaxy formation and evolu- 
tion from our own ‘backyard’, the Milky Way galaxy and its satel- 
lite system. Resolved stellar spectroscopy of the Milky Way stellar 
halo and the satellite galaxies provides ‘archaeological’ evidence 
of the chemical enrichment of the interstellar medium back to the 
earliest times. The details with which the Milky Way and its satel- 



lites can be studied make them a useful testbed of the cosmological 
paradigm. 

The star formation history of the Milky Way satellites can be 
derived from the study of stellar populations identified in co lour 
magnitude diagrams (CMD) (see iTolstov. Hill & Tosill2009L and 
references therein). These studies have revealed large variations in 
the star formation histories of Local Group dwarf galaxies, even 
for those of similar stellar mass. These range from solely old- to 
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predominantly intermediate-age to even significantly young stel- 
lar populations. Some st ar formation histories may be bursty, s uch 
as for the Carina dSph iHurley-Keller, Mateo & Nemedll998ll . It 
is currently not completely understood what physical mechanisms 
are responsible for the exact star formation history and how the en- 
vironment of the galaxy influences the star formation process. 

With additional spectroscopic observations, several teams 
have investig ated the dynamical and chemical proper t ies of both 
classical (e.g. TolstoveLal 20061: Battaglia et^ l2006l : lKoch et ^ 



l2007bl : IWalker et alJl2009F Kirbv et alj|20ici|) ' and the recently dis- 
covered ultra-faint dwarf spheroidal gala xies (e.g., IKirbv et aO 
l2008l : IXden et aljr2009l : iNorris et al. 2010l). The discovery of this 



new, v ery faint class of satellites ( WlhnMiet^ .120051:1 |Zuckeret^J 
2006a|[^; jBek^m et al.l l2007l : IWalsh. Jerien & WillmanI l2007l : 



Irwin et al ■11200/1) h as revived the interest in the so-called ‘ missing 



satellites’ problem iKlvpin et alJI 19991 : iMoore et alJll999l) . which 
contrasts the huge number of dark matter satellites predicted to or- 
bit in Milky Way-sized haloes with the relatively modest number 
of luminous satellites observed. 

In parallel to observational efforts, large cosmological N- 
body dark matter on ly simulations, l ike the Aquarius Pr oject 
dSpringel et al.ll2008a|). Via Lactea II IPiemand et alj IToOSl) and 



GHALO I Stadel et^l2009l) . have greatly improved mass resolu- 



tion and have now reached a regime in which the formation and 
evolution of (satellite) galaxies can be studied in exquisite detail 

down to this ultra-faint regime in a ACDM universe. 

Following the early sugges ti on by lEfstathio^ lll992h and 
Kauffmann, White & Guiderdonil ( 1 19931) that the reionization 
of the intergalactic mediu m at high redshift could suppr e ss the 
formation of faint galaxies, Bullock, Kravtsov & Weinberg] j2000t) 
were able to show, using dark matter halo merger trees and a 
1-D gas simulation, that the effects of reionization could indeed 
help to reconcile the distribution of subhalo circular velocities 
expec ted in the CDM co smology with inferences from satellite 
data. iBenson et al.l j2002h then developed a detailed treatment 
of reionization and, using a semi-analytic model of galaxy for- 
mation, showed that the combined effects of reionization and 
supernova feedback could account for the observed luminosity 
function of satellites in the Local Group. They also predicted 
the existence of a large population of ultra-faint satellites. Sev- 
eral recent semi-analytical and hydrodynamical studies have 
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back mechanisms, to suppress the formation of small galaxies 
within all haloes and reproduce the observed number of dwarf 
satellites down to the ultra-faint regime. 

However, many issues remain. Relevant questions to be asked 
are for instance: ‘How many satellite galaxies are still undis- 
cover ed in the Milky Way stellar halo? ’ (e.g. IKoposov et al' 



2008 ), ‘What was their time of infall?’ (e.g.. lLi, De Lucia & Helmi 



‘What is the mass of 



the Milky Wav dark matter halo?’ (e.g 


Battaglia et al.l 20051: 


Sales et alj 2007 


: Smith etalJl200llLi& White 2008: Xue et all 


20081: IGuo et al 
lites linked to 


I 2 OIOI). ‘Are the luminosity functions of satel- 
he properties of their host in any way?’ (e.g. 


McConnachie & IrwinI 20061: McConnachie et alj 20091: Guo et alj 


201 lal:lLares. Lambas & Dommguezl201 iMWang & Whitel2012al). 



In this work, we study the formation and evolution of dwarf 



galaxies in and around Milky Way-lik e galaxies using the N -body 
simulations of the Aquarius Project ISpringel et al. l2008ah . We 
combine these with semi-analytical modelling to study the physical 
processes associated with the bary onic components of the galax - 
ies. We use the model described bv iLi. De Lucia & Helmil ( l2010h . 
which has been extended to include new prescriptions to follow the 
stell ar stripping and tid al disruption of satellites. 

iFont et ai] (1201 ih also combined the Aquarius simulations 
with a semi-analytical code to study the properties of satellite 
galaxies, but their focus was in particular on a more sophisticated 
treatment of reionization, while our interest is mainly in the star for- 
mation histories of the satellites and of isolated dwarf galaxies. Be- 
cause the two codes were developed independently, it is instructive 
to compare their results on general properties for the satellite pop- 
ulation, such as luminosity function, metallicity distri bution and 
radial profile. A sim ilar semi-analytical m odel to that of iFont et alj 
l201lh was used by ICooper et al.l ll2010l) to study the stripping of 
satellite galaxies and the format ion of Gala c tic halo es in the Aquar- 
ius simulations. Additionally, I Guo et all d201ld) have used an 
adapte d version of the semi-analytical code of lPe Lucia & Blaizod 
d2007h to study galaxies and satellites in the Millennium II simula- 
tion and show their results to be consistent with the satellite lumi- 
nosity function over the (lower) resolution range in that simulation. 

This paper is structured as follows. In Sections 12.11 and 12.21 
we describe the Aquarius simulations and the particular semi- 
analytical model we use. Some additional prescriptions have been 
implemented to account for the tidal stripping and disruption of 
satellites, and these are described in detail in Appendix lAl In Sec- 
tion we investigate the properties of the modelled main galax- 
ies as well as the luminosity function, the luminosity-metallicity 
relation, and the radial and spatial distributions of their satellites. 
Section |4] is devoted to a more in-depth analysis of the star forma- 
tion histories of the modelled dwarf galaxies, both satellites and 
isolated galaxies, whereas in Section |5] we investigate the clos- 
est model analogs to the dwarf spheroidal galaxies Sculptor, Ca- 
rina and Fornax. In Section we discuss our findings and com- 
pare them with other semi-analytical and hydrodynamical work. 
We summarise our results in Section|7] 



2 THE MODEL 

2.1 The Aquarius simulations 



The six Milky Way-like haloes (Aq-A to Aq-F) of the Aquar- 
ius Project were selected from a lower resolution ver sion of the 
Millennium-II Simulation iBovlan-Kolchin et alj|200^ ). a cosmo- 
logical N-body simulation of a cubic region 125h“^Mpc on a side 
with parameters Dm = 0.25, Ha = 0.75, ag = 0.9, ris = L 
h = 0.73 and Ho = lOO/i km s“^Mpc“^. We refer the reader 
to lSpringel et al.l ll2008al l9) for further information. The parameters 
are th e same as those of the Millennium Simulation dSpringel et alj 
l2005h and were based on the first-year results from the WMAP 
satell ite. They are no longe r consistent with the latest WMAP anal- 
ysis iKomatsu et alj l201lh . but we do not ex pect this to affect 
our results significantly (see IWang et al.lE008L for a comparison 
of first and third year parameters). The simulated Milky Way- like 
haloes have virial masses (M 200 , defined as the mass enclosed in 
a sphere with mean density 200 times the critical value) in the 
range 0.8-1.9X lO^^M©, broad ly consistent with the mass esti- 
mated for the Mill ^ Way (e.g. Battaglia et al. 20051: ISmith et all 
I2OO7I: iLi & Whitel I2OO8I: IXue et alJ l2008l: IGuo et all 120101) . One 
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Name 


mp 


M 200 


^200 


Nr. snapshots 




(M©) 


(M©) 


(kpc) 




Aq-A-5 


3.14 X 10® 


1.85 X 10®® 


246.37 


128 


Aq-A-4 


3.93 X 10® 


1.84 X 10®® 


245.70 


1024 


Aq-A-3 


4.91 X 10 “ 


1.84 X 10®® 


245.64 


512 


Aq-A-2 


1.37 X 10 “ 


1.84 X 10®® 


245.88 


1024 


Aq-B-2 


6.45 X 10® 


8.19 X 10®® 


187.70 


128 


Aq-C-2 


1.40 X 10 “ 


1.77 X 10®® 


242.82 


128 


Aq-D-2 


1.40 X 10 “ 


1.77 X 10®® 


242.85 


128 


Aq-E-2 


9.59 X 10® 


1.19 X 10®® 


212.28 


128 


Aq-F-2 


6.78 X 10® 


1.14 X 10®® 


209.21 


112 



Table 1 . Some basic parameters for the Aquarius haloes from ISpringel et all 
l2Q08ah (see the original paper for more information). The columns corre- 
spond to the simulation name, the particle mass (rrip), the virial mass of 
the halo (M 200 ). Ihe corresponding virial radius (r2oo) the number of 
snapshots we use for each simulation. 



halo, Aq-A, was simulated at five different numerical resolution 
levels (summarised in Table [T). We focus on the high-resolution 
level 2 (common to all six haloes) and use the lower resolution se- 
ries Aq-5 through Aq-2 to test the numerical convergence of our 
model. 

Dark matter haloes are identified in the simulations using a 
friend s-of-friends algorithm jPavis et al.lll985h and the code SUB- 
FIND dSpringel et aPl20(^ which identifies self-bound structures 
within larger structures. Following previous work, we have only 
considered subhaloes that retain at least 20 particles. 



2.2 The semi-analytical code 



We use the Aquarius simulations as a backbone for modelling bary- 
onic processes in galaxies. Subhalo catalogues are used to construct 
merger (history) trees for all self-bound haloes and subhaloes in 
the A quarius simulations llSpringel et al.l200^lDe Lucia & Blaizoll 
l2007h by determining one unique descendant for each (sub)halo. 
These merger trees are combined with semi-analytical modelling 
to study the galaxies that reside in such subhaloes. The semi- 
analytical modelling technique follows the relevant physical pro- 
cesses using simple but observationally and astrophysically moti- 
vated ‘prescriptions’. One advantage of this method is that it is ap- 
plicable to large cosmological simulations and provides relatively 
fast predictions of galaxy properties. The method, however, does 
not follow explicitly the gas dynamics (as is done in hydrodynam- 
ical simulations), and does not usually provide spatially resolved 
information about the baryonic components. 

The sp ecific model we use in this wo rk is the ‘ejection model’ 
described in lLi, De Lucia & Helmil (l2010h . with new prescriptions 
to follow the stellar stripping and tidal disruption of dwarf galax- 
ies as they become satellites. These new prescriptions are de- 
scribed in Appen dix IaI The mode l build s upon the me t hodol- 
ogy introduced by |K aufflnann_gLffl \ \9% ), ISpringel et alj j200lh 
and De Lu cia. KauffiMnn_& White (2004 ) and subseque n tly up - 
dated by ICroton et al.l (2 0061) and De Lucia & Blaizol 



The model of |Pe Lucia & Blaiz^ j2007h has since been mod- 
ified to follow more accurately p rocesses on the s c ale o f 
the Mil ky Wa y and its satellites by D e Lucia & Flelmil (20081 ): 
Li et alj 1 2009h: iLi, De Lucia & HeliS ( 2 OI 0 I) . In particular, both 



Li, De Lucia & He lmi|(201 Cjl and this paper use the ‘ejec tion’ feed- 



back scheme of |Pe Lucia. Kauffmann & Whitd (20041 ). which is 



different from the (default) feedback s c heme adopted in most pre- 
vious work (e.g., icrotqn et al.l I2OO6I : |Pe Lucia & Blaiz^ 120071 
|Pe Lucia & Helmil |2008|) . a somewhat earlier reionization epoch 
and additiona l ly supp ress cooling in small h aloes (Tyir < 10^ 
K). Guo et all (201 la ) show that the model of lPe Lucia & BlaizoJ 
(2007!) significantly overpredicts the number of galaxies with stel- 
lar masses between 10^ and IO^^Mq. Our model is based on 
|Pe Lucia & Blaiz^ (20071) . but as outlined above several changes 
have been made since that paper. Most changes are focused on im- 
proving the modelling on low mass scales and are shown to mostly 
affect the dwarf galaxy scale and preserve th e properties of galaxies 
with a Milky Way-like mass (see Table T_of[Lij^£Luci^&HeM 
I2OI0L for a comparison). However, |Pe Lucia & Blaizol ~(2007l) 
show that the choice of feedback scheme will also affect the lu- 
minosity and evolutionary rate of t he brightest cluster galaxies. In 
particu lar, the feedback model of |Pe Lucia. Kauffmann & White! 
(2OO4I1 . as used in this work, will result in a more prolonged star 
formation activity and a higher luminosity f or these mass i ve gal ax- 
ies compared to the feedback model used in ICroton et al.l (20061) . A 
careful analysis of the combined impact of all the model changes 
on more massive galaxies will be part of future work. 

A schematic diagram of the main processes modelled, is 
shown in Figure [T] Below, we outline the main ingredients of the 
model, but for a fu ll description and the analytic expressions we 
refer the reader to iLi, De Lucia & Helmil ( 2 OIOI) and references 
therein. 



• Reionization is modelled following iGnedinI (2000l) . and 
ICroton et al.l (20061) . Reionization causes the baryonic content of a 
halo to decrease in haloes with mass comparable or smaller than a 
so-called filtering mass, which evolves with redshift. The reioniza- 
tion epoch is assumed to last from zq = 15 to = 11.5. Our reion- 
ization prescription result s in a stronger effect than sugge sted for 
the g lobal reion i zation in lOkamoto, Gao & TheunsI (20081) . How- 
ever, iFont et aP (201 ll) show using a detailed treatment of reion- 
ization for the Aquarius simulations that the proto-Galactic region 
is completely photo-ionized by 2 = 10 due to the contribution of 
local sources. In the end, this results in a comparable reionization 
history on t he scale of the M ilky Way to that in our model (see Ap- 
pendix A of iFont et alJl201 iL and references therein for a complete 
discussion). 

• Cooling of the hot gas is dependent on its metallicity and tem- 
perature. Below Tvir = 10^ K (the atomic hydrogen cooling limit) 
cooling is forbidden since we are assuming that cooling via molec- 
ular hydrogen is prevented by photo-dissociation caused by UV ra- 
diation from the (first) stars in most cases (in contrast to other semi- 
analytical work where star formation in haloes with Tvir < 10'* K 
is allowed through molecular hydrogen cooling; see for instance, 
iMadau et alj|200a ; ISalvadori, Ferrara & SchneideilbOO^ . 

• Star Formation transforms cold gas, which is assumed to be 
in an exponential thin disc w ith properties given by the formalism 
of lMo, Mao & Whitel l fl^ . into stars. Stars form in the gas of the 
disc that is above a critical density threshold. The star forming disc 
radius, idisc, is assumed to be three scale-lengths. Assuming that 
the disc has a flat rotation curve with a rotation al velocity equal to 
the circular velocity of the halo (V200, see also iKauffmannll 199a) 
and the gas velocity dispersion is 6 km s~^, the cr itical density 
threshold is described by Eauation[T1 (Kennicutlll 19891) . 



^crit „ rrv ^200 / fdisc 

^ IT = L).5y ^ 

M0PC-2 kpc 



(1) 



The star formation rate is then proportional to the amount of gas in 
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Figure 1. A schematic diagram of the semi-analytical modelling scheme for the central galaxy within a dark matter halo and a satellite that has just been 
accreted. The yellow boxes linked to the galaxies represent all different ‘phases’ of the baryons (these are all modelled analytically) and the red arrows 
represent all modelled physical prescriptions that affect them. 



this state. Star formation can also happen in bursts during minor or 
major mergers, when (part of) the cold gas of the mergi ng galax- 
ies is turned into stars. A Chabrier initial mass function (IChabrien 
l2003h and instantaneous recycling approximation are assumed and 
accordingly 43% of the mass in stars formed at each time step is 
(instantaneously) recycled back into the gas phase, representing su- 
pemovae Type II events. 

• Heating and ejection of gas is due to supernova feedback pro- 
cesses. We u se the feedback prescription from the ejection model 
described in |Pe Lucia, Kauffmann & Whitd | |2004| ). The mass of 
gas reheated by supemovae depends on the depth of the halo poten- 
tial well (i.e., (X l/V 2 oo)i which implies that smaller haloes (with 
shallow potential wells) are more sensitive to the effects of feed- 
back. The material reheated is put in an ejected component that can 
be reincorporated into the hot gas at later times. 

• Metals are created in star formation events and follow the flow 
of the mass between different components. We assume an instanta- 
neous recycling approximation that is appropriate only for elements 
produced by Type II supemovae. 

• If a galaxy becomes a satellite (i.e., its halo becomes a sub- 
halo) we assume that its hot gas and ejected gas components 
are stripped and transferred to the central galaxy. This process 
cmdely models the physical process of ram-pressure stripping of 
hot gaseous components. In this paper we have also added two 
physical mechanisms that only operate on satellite galaxies due 
to interactions with the host halo, stellar stripping and tidal dis- 
ruption. The implementation of these physical processes and their 
results are described and shown in Appendix We find that our 
implementation of stellar stripping does not have a significant ef- 
fect on the luminosity functions, it affects very few galaxies. The 
tidal dismption prescription allows us to deal with galaxies which 
have ‘lost’ their dark matter subhalo in the simulation (when it is 
stripped down to fewer than the SUBFIND resolution of 20 particles) 
and decide whether the galaxies themselves should survive the tides 
of the main galaxy. We refer to these galaxies as ‘orphans’, in con- 
trast to satellites which still live in dark matter subhaloes of more 
than 20 particles. Our implementation of tidal dismption has a sig- 
nificant influence on the shape of the satellite luminosity function, 
as we will show in Section lT^ 



3 COMPARISON TO THE MILKY WAY AND ITS 
SATELLITES: GENERAL PROPERTIES 

3.1 Milky Way properties 

In Figure 12 we show the various properties of the central (Milky 
Way-like) galaxies residing in the main Aquarius haloes compared 
to the values ob served for the Milky Way. W e show the results of 
our basic model JlI, De Lucia & HelmM201(X filled circles) and the 
same model with our additional prescriptions for stellar stripping 
and tidal dismption (see Appendix lAl open squares). Our treatment 
of satellite stripping and dismption changes only slightly the prop- 
erties of the main galaxies. The scatter from main galaxy to main 
galaxy in the different Aquarius haloes is clearly visible. 

Aquarius galaxies B and E have very similar stellar masses 
to the Mi lky Way, estimate d to be in the range of 4.9 — 5.5 x 
10^° M0 jFlvnn et alj|2006l) . one of the most robust and important 
constraints for which the comparison can be made. 

In our models the spheroid refers to both the bulge and stel- 
lar haloes of the galaxies. Since in the case of the Milky Way 
the stellar halo contains a very small mass (in comparison to the 
bulge), we compare the spheroidal component in the model with 
the bu lge component of the Milk y Way. As explained in more de- 
tail bv IPe Lucia & Helmil j2008l) . spheroid formation can occur in 
our model through both mergers and disc instabilities. The treat- 
ment of disc instabilities is one of the least constrained physical 
processes in the model. It is very sensitive to small variations of the 
other prescriptions, but nevertheless is a key channel for bulge for- 
mation in Milky Wa y size galaxies (see lParrv, Eke & FrenklboO^ : 
|Pe Lucia et alj2oIIL for a comprehensive discussion). The disc and 
bulge components of the M ilky Way are estimated to have a mass 
ratio between 0.2 and 0.3 iBissantz, Pebattista & GerhardI l2004l) 
which implies a bulge mass of 0.8 — 1.3 x 10^° Mq. This bulge- 
to-disc ratio is well reproduced by the Aquarius galaxies, except 
in the case of galaxies E and E, which have a very dominant bulge 
component. 

The slight increase of the total spheroid mass (and thus also 
of the bulge/disc ratio) in the models with stripping and tidal dis- 
ruption is expected since the stars stripped from the satellites are 
added to this component. Additionally, these processes can affect 
the details of the disc instability phenomenon through an increase 
in the cold gas mass of the disc. However, we find that for all haloes 
the major source of the bulge mass increase is the addition of stel- 
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Figure 2. Properties of the central ga laxies of the haloes Aq-A to F both 
for the ‘satellite-model ejection’ from iLi. De Lucia & Helm] i2QlQh (filled 
circles) and our fiducial model (i.e., extended to include stellar stripping and 
tidal disruption, open squares). The solid lines show observed Milky Way 
values (see text for references). The properties of the Milky Way bulge are 
compared to the spheroid component in our models. 



lar mass from disrupted and stripped satellites. For instance, the 
relatively large difference found in Aq-F is almost completely due 
to one quite lumin ous satellite galaxy that is tidally disrupted (see 
ICooper et al.ll201 iL which includes a movie of the evolution of this 
object). 

In all models, except Aq-E which shows almost no star forma- 
tion at the present day, we find that the current star formation rate is 
above the corresponding value for the Milky Way. This is also true 
for the total amount of present-day cold gas, which is also above 
the Milky Way value for most models (with the exception of Aq-F) 
(using a total HI-)-H 2 mass value for the Milky Way of ~ 6 x 10® 
Mq, see lBlitjll997l and references therein). Among the models 
showing significant ongoing star formation, halo B and F give the 
closest match to the Milky Way. 

The two bottom panels of Figurel^show the average metallic- 
ities obtained for all the stars and just the spheroid component. The 
overall metallicities lie systematically above the Milky Way values, 
but despite the crude approximations made to follow the evolution 
of metals in our model, the mismatch is only on average ~0.25 
dex, or a factor ~1.8 in total mass of metals. Galaxy B, E, and F 
show the best match to the measured [Fe/H] for the disc (0. 1-0.2 
dex above the measured value), and halo B a very close match to 
the bulge stars jFreeman & Bland-Hawthornll2002h . 




-5 -10 -15 -20 

Mv 



Figure 3. Cumulative luminosity functions of all satellite galaxies within 
280 kpc of the main galaxy for all different Aquarius haloes (solid lines), for 
the Milky Way satellites (black filled circles), with corresponding Poisson 
error bars, and for the Milky Way satel lites as derived and corrected for 
incompleteness by IKodosov et al J l2008li (thick dashed black line). In the 
top panel surviving orphan satellite galaxies are also included, in the bottom 
panel they are not shown. 



Altogether, Figure [^suggests that the galaxy modelled within 
halo B is the closest analogue of the Milky Way galaxy. 



3.2 Satellite luminosity functions 

Figure[3]shows the luminosity function of all satellites of the main 
Aquarius galaxies (solid lines of different colours). Black circles 
show the cumulative number of Milky Way satellites. Given that the 
Sloan Digital Sky Survey (SDSS) footprint covers roughly a quarter 
of the sky, we may expect many more galaxies hiding in the Milky 
Way that are either outside the SDSS footprint or t oo faint and/or 
too fa r away to be identified in SDSS data. The IKoposov et alj 
d2008l) (thick dashed black line) relation takes these uncertainties 
into account and attempts to correct for them. However, it should be 
noted that for the brightest end (Mv < — ID IKodosov et alld2008l) 
use an average luminosity function of the Milky Way and M3 1 for 
their fit. Since M31 has more bright satellites than the Milky Way, 
the relation overpredicts satellites in this regime. 

In the top panel of Figure[3]both the satellites still embedded 
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within a dark matter subhalo and the orphan satellites (for which 
their dark matter subhalo has been stripped below the resolution 
limit of the simulation) which survive according to our tidal disrup- 
tion prescription (see Appendix lAt are plotted. In the bottom panel 
these orphan satellites are not shown. Both the data and the model 
counts of satellite galaxies are restricted to a dista nce of 280 kpc 
from t he centre of the main galaxy, as assumed by IKoposov et alj 
d2008h . This radius is comparable to, but generally a bit larger than, 
the virial radii (r 2 oo) of the main Aquarius haloes (see TablelTJ. 

The shape of the luminosity function roughly agrees with the 
Milky Way data down to Mv = — 5 where resolution effects are 
starting to play a role (see Appendix m for the discussion of this 
limit, the vertical dashed line in Figure |3 for all Aquarius haloes. 
Halo B shows clearly the best quantitative correspon dence with the 
luminosity function derived bv lKoposov et al.l d2008l) in our fiducial 
model (top panel). It is interesting that Aquarius B has a similar 
satellite luminosity function as well as the central galaxy that most 
closely resembles the Milky Way galaxy, as shown in the previous 
subsection. However, one should bear in mind that the number of 
satellites formed around any halo is very sensiti ve to the choice of 
feedback schem e and reionization physics (e.g.. lGuo et alJl^I Ibl : 
iFont et alj20Ilh . 



Most Aquarius haloes do not co ntain satellites as brigh t as the 
LM C. This has been i nvestigated by iBovIan-Kolchin et alJ d2010l) 
and iBusha et ^ J2^ , who estimate from dark matter simulations 
that the probability of finding both an LMC and SMC around a 
Milky Way-sized halo is ~ 10 per cent (up to 25 per cent with a 
dependence on the exact Mil ky Way dark matter halo mass and en- 
vironment). fkiu et al.l d201 ih find from an analysis of Milky-Way- 
like hosts in the SDSS DR7 catalogue that only 3.5 per cent of them 
have two such bright satellites within 150 kpc of their host. A fur- 
ther warning against making too strong a st atement about the bright 
end of the luminosity function comes from iGuo et al.l d201 lah . who 
find that isolated host galaxies of luminosity comparable to the 
Milky Way and to M3 1 c ontain ^2 times fewer satellites brig hter 
than My = — 14 (see also lLares, Lambas & Dommguezll^Ol ih . In 
our models, 4 out of 6 haloes have one satellite as luminous as or 
slightly more luminous than the SMC (the exceptions are haloes B 
and C). Aq-F hosts a satellite galaxy slightly brighter than the LMC 
and one of luminosity similar to the SMC. Halo D hosts three very 
luminous satellites, but all fainter than the LMC. 

Although all Aquarius haloes have masses consistent with that 
of the Milky Way halo, they still span a factor 2.25 in mass. Aquar- 
ius B is the least massive with M 200 = 8 x 10^^ M©. Figure |3 
shows that this range of mas ses is reflected in th e number of satel- 
lites. This was noted also bv iMaccio et alj ( I2OIOI1 . who remark that 
the trend between halo mass and satellite luminosity function does 
not depend on the particular semi-analytical model used. For our 
semi-analytical model this correlation is shown in Figure [4] where 
the number of satellites brighter than Mv=-5 within the virial ra- 
dius is plotted against the virial mass of the host halo. The trend 
is very clear and the scatter is small, although based on just a few 
points. We plan to investigate these results further using the much 
larger volume of the Millennium II simulation, where over 70 00 
Milky Way-sized haloes are found jBovlan-Kolchin et alj|201(il) . 

Since our model for halo B best resembles both the properties 
of the Milky Way and the luminosity function of satellite galaxies, 
this would favour a dark matter mass estimate close to 8 xlO'^^ 
M(7) for the Milky Way galaxy in agreement with the wor k of 
iBattaghu eUiyj ^ : ISrnith et^ j2()07ll : I Sales et ^ l l2007h and 
Ixue et alTi 2008|), bu t a fac tor t wo less massive t han the best esti- 
mate of iLi & Whitel l l2008h and iGuo et all j201(]|) . In terms of for- 




Figure 4. The virial masses of the six different Aquarius dark matter haloes 
versus the number of satellite galaxies of the main galaxy in our model. 
Open black squares show the number of satellites brighter than Mv=-5 
within the corresponding virial radius. 




Mv 

Figure 5. Luminosity and metallicity for the satellite galaxies (grey filled 
circles) and those present only in halo B (black filled circles). Grey open 
circles represent satellite metallicities we do not trust, mainly due to in- 
complete modelling of first star physics. Overplotted as red asterisks are 
the average values for the Milky Way satellites, corrected to approximate 
a mass-weighted average of [Mg/H] for a better comparison to the models 
(see text for details and references). The error bars indicate the metallicity 
scatter found inside the galaxies. 



mation history, halo B forms relatively late, its mass accretion is 
slower than that of other Aq uarius haloes, in particular at 2 : > 2 
teovlan-Kolchin et alj|20ic|) . In other properties, like spin o r con- 
centration, halo B is not special teovlan-Kolchin et alj201o|) . 
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3.3 Luminosity-metallicity relation 

The Milky Way satellites (including the ultra-faint dwarf 
spheroidals: iKirbv et alll200^ show a strong correlation between 
luminosity and metallicity. It is not straightforward to compare 
the observed abundances to the metallicities in our model, since 
the model adopts an instantaneous recycling approximation that is 
valid for the majority of alpha elements formed like O and Mg, 
but not for Fe (mainly produced by Type la SN) which is the most 
commonly measured element in stellar spectra. 

For our comparison we therefore use Mg. Correcting the ob- 
served average [Fe/H] to [Mg/H] requires knowledge of [Mg/Fe]. 
Observations of red giant branch stars in dwarf galaxies and the 
Milky Way stellar halo show clear trends of [Mg/Fe] with [Fe/FI], 
which vary slightly from satellite to satellite and are distinct from 
those in the Milky Way stellar halo, especially at [Fe/H]> —1.5. 
Nonetheless, we adopt th e following function based on observed 
satellite data compiled bv lTolstov. Hill & Tosil j2009h for [Mg/Fe] 
with [Fe/H]: 



their few star formation episodes these satellites will enrich very 
little and not compensate for the first generation of stars formed 
with no metals at all, in contrast to higher mass objects which can 
sustain star formation for more extended periods. Their very low 
metallicities could therefore result mainly from our neglect of any 
kind of pre-enrichment, which could be driven by a top-heavy IMF 
for the first stars. These events are likely to enrich the galaxy, or 
even the intergalactic medium to a metallic ity floor of [Fe/H] ^ — 3 
(e.g. lSalvadori, Ferrara & SchneiderboOSh . 

The metallicities for SMC and LMC-like satellites are over- 
predicted, although they are found to be consistent within the range 
of the measured scatter inside these galaxies. A more thorough 
modeling of the chemical processes will have to be conducted how- 
ever, before enabling any conclusions on the possibility that our 
model retains too much of its produced metals for the more massive 
galaxies. The overall slope and normalisation of the metallicity- 
luminosity relation are reproduced well in our models when con- 
sidering those galaxies which have sufficient star formation events 
to enrich above [Fe/H]=-3, as shown in the grey filled circles. 



[Mg/Fe] = -1-0.4 
[Mg/Fe] = -0.4[Fe/H] - 0.4 
[Mg/Fe] = -0.4 



for [Fe/H] < -2 
for - 2 < [Fe/H] < 0 
for [Fe/H] > 0 



(2) 



The metallicity given by our model is mass-weighted, i.e. it 
is the logarithm of the ratio of mass in metals over the total mass 
in stars. The average observed metallicity to compare with should 
therefore also be obtained by taking the logarithm of an average 
over the ratio of metals to hydrogen, which will be different from 
the average of [Fe/H]. We have tested the offset between the two es- 
timates of the mean metallicity on the data se t of the Dwarf Abun - 
dances and Radial velocity Team (DART) iTolstov et alJ l2004h . 
which contains (Ca II triplet derived) metallicities for the classical 
dwarf spheroidals Fornax, Sculptor, Carina and Sextans. We found 
that the average [Fe/H] is on average 0.23 dex lower than the log- 
arithm of the average over the ratio of Fe to H, although ranging 
from 0.15 to 0.4 dex. In Figure we therefore show [Fe/H] cor- 
rected by 0.23 dex and subsequently transformed to [Mg/H] values. 
For Fornax, Sculptor, Sextans and Carina we do not use the 0.23 
correction, but use t h e averages dire c tly fro m th e DART data (see 
iTolstoy et al.l ( l2006h . [Battaglia et alj j2008l) and IStarkenburg et alj 
J 2 OI 0 I) for a description of the data sets and methods) by taking the 
logarithm of the average ratio of Fe to H for all stars. 

Figure |5] shows a comparison for those Milky Way satel- 
lites which have an average [Fe/H] avail able f rom the literature 
(see, for instance, the compilation in Ol200^ . The mean iron 
abun dances and their di spersions (error bars) in t he plo t are ta ken 
from lWesterlundl lll997h for t he LMC and SMC, IColel d200lh for 
Sagit tarius, the DART survey llHelmi et alj2006l : lstarkenburg et alj 
|201(1 and reference s there in) for Fornax, Sculptor, Cari na and Sex- 



tans. 



Harbeck et alJ ||2001 ) for Draco and Ursa Minor, &chet^3 



for Leo I and Koch et aT] l l2007j) for Leo II, iKirbv et^ 



1 20080 for most of the ultra-faints; Bo otes I and Segue I are from 
iNorris et al.l ( 1201 Ot) and Bootes II from iKoch et alj d2009h . 

The filled and empty grey circles in Figure|5]show the average 
metallicity for the model satellite galaxies within 280 kpc of their 
hosts. The larger black filled circles represent the satellites in Aq- 
B, to highlight the number of satellites and dispersion in metallicity 
found in that simulation (which most resembles the Milky Way). 

All model galaxies with an average metallicity below [Fe/H] 
=— 3 are shown as open grey circles. All of them have experienced 
very few star formation events (typically less than 4). Because of 



3.4 Radial and spatial distribution 

In Figurel^we show the radial distributions of all satellites brighter 
than Mv=— 8.5 within the Milky Way and in each of the Aquar- 
ius haloes. For the (surviving) orphan satellites included in the top 
panel, the position is that of the most bound particle within the host 
subhalo before disrup tion. For the M ilky Way satellites, the dis- 
tances are taken fr om iMateol dl998ll . except for Canes Venatici I 
dMartin et alJl200^ . We have assumed a distance from the Sun to 
the Galactic centre of 8.5 kpc. The uncertainty in the radial distri- 
bution due to Poisson noise is indicated by the grey area. 

The distribution of satellites is broadly consistent with that ob- 
served in the Milky Way for most Aquarius haloes, except for the 
inner regions. If orphan satellites are removed, the profiles are in 
general slightly less centrally concentrated, since the orphan galax- 
ies are mainly found in the inner regions. The bottom panel of Fig- 
ure shows that the radial distributions of Aquarius C, D, E, and 
F match that of the Milky Way if orphan galaxies are excluded 
from the analysis. Such an exclusion could be justified, since or- 
phan galaxies are almost certainly affected by tides and could have 
very low surface brightness, which would hinder their detectability. 

Note that the radial distribution of luminous satellites is dif- 
ferent from that of the total population of dark matter substruc- 
tures, which has a less centrally concentrated pro file as shown as 
dashed lines in the to p panel of Figure (see also lGao et al.ll2003 ; 
ISpringel^ 'Zll2008ah . 



One interesting property of the Milky Way satellites is 
that they seem to lie close to a plan e, rather than being dis- 
tributed isotropica ll y on the sky (e .g. iKunkel ^.Dgrngri [l97^; 



Lvnden-Bell 19761; 


Maiewskil 1994 Lvnden-Bell & Lvnden-BeUl 


I 995 I; HartwickI 


200d; IPalma. Maiewski & Johnstoij |2002|; 


Krouoa. Theis & Boilvl 20051; [Metz. Krouna & JerienI 2007h. In 



Figure |7] we investigate the anisotropy of the spatial distribution 
of satellites with Mv < —8.5 in each of the Aquarius haloes. We 
calculate the flattening using the normalised inertia tensor. The 
short-to-long axis ratio (c/a) is computed from the eigenvalues 
of the diagonalized inertia tensor. In each panel the flattening of 
the Milky Way satellite system is indicated by a dashed black 
vertical line, the flattening of the total system of subhaloes in each 
Aquarius halo as a dashed grey line, whereas that of all bright 
satellites is shown as a full coloured line. 

All Aquarius haloes show a fairly spherical distribution of 



© 0000 RAS, MNRAS 000. 000-000 



8 E. Starkenburg et al. 




R (kpc) 

Figure 6. Radial distribution of the bright model satellites (My < —8.5, 
solid lines). The top panel shows all bright satellites (solid lines) as well as 
the total population of all subhaloes (dashed lines). The grey area indicates 
the distribution of Milky Way satellites including a Poissonian error bar 
in both panels. In the bottom panel additionally the radial distribution of 
classical Milky Way satellites is overplotted (dotted black line) and orphan 
satellites are excluded from the solid lines. 

their total system of (dark) subhaloes within 280 kpc, although 
some variations can be seen from halo to halo. The spatial distribu- 
tion of their bright satellites is in all Aquarius haloes less flattened 
than the Milky Way satellite system. However, all host a larger 
number of bright satellites as well, as shown in Figure |3 To in- 
vestigate the effect of the number of satellites on this comparison, 
we have overplotted in all panels in Figure |7] the distribution of 
c/a when instead of all bright satellites a random subsample of 12 
satellites is taken, equal to the number of classical satellites in the 
Milky Way. The restriction to a smaller number of satellites greatly 
enhances the chances of selecting a more flattened distribution, as 
can also be concluded from the bottom panel where 12 points are 
selected randomly distributed on a sphere. From a purely spheri- 
cal distribution one expects a flattening comparable to that seen in 
the Milky Way in ~1 per cent of the cases if 12 points are drawn 
randomly. 

The chance of getting such a highly flattened distribution as 




c/a 



Figure 7. Flattening (c/a) of the twelve classical satellites (My < -8.5) in 
the Milky Way (black dashed vertical line) and the distiibution of c/a values 
that can be reached with random selections of twelve My < -8.5 model 
satellites for the Aquaiius haloes A-F. For each Aquarius halo also the flat- 
tening of all satellites with My < -8.5 is overplotted as a full coloured line. 
We also show the flattening of all subhaloes within 280 kpc of each main 
halo as a dashed grey line. The distribution of c/a obtained by taking twelve 
random points from a spherical distribution is shown in the bottom panel 
for comparison. 



seen in the Milky Way from 12 satellites within the Aquarius mod- 
els is low, but can also not be completely ruled out. In particular, 
Aquarius B and E show distributions close to the spherical case 
when 12 satellites only are selected. It is extremely unlikely for 
Aquarius B to host a Milky Way-like flattened system (but notice 
that Aq-B has significantly fewer bright satellites in total than the 
other haloes). Aquarius A, C and F are on average more flattened. 
In some of these cases the flattening of the satellite system fol- 
lows the shape of the present-day host dark matter halo and/or large 
scale structure, most notably in Aquarius A as illustrated in Figure 
[8] Aquarius A is found to have a long, thin filament which is co- 
herent in time, whereas in some of the other Aquarius haloes the 
filament is either less well-defined or broader (such that it encom- 
passes the whole halo) or changes its orientation over time (for a 
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Figure 8. The spatial distribution of galaxies of —19 <Mv < —8.5 in 
Aquarius A. Symbols and colours denote the age of the dominant stellar 
population: old (red asterisks), or intermediate (blue diamonds). The grey 
map shows the underlying density of dark matter in the same simulation. 
The frame is rotated such that the major axis of the main halo is vertical 
and the minor axis horizontal. The grey errors show the relative sizes of the 
major and minor axis. 



full discussion of the shapes of the Aq uarius dark matter ha loes 
and their filaments we refer the reader to lVera-Ciro et al. 11201 ih . 

Other studies have also indicated that a flattening similar 
to that of the Milky Way satellites can be reach ed in a ACDM 
cosm o logy, although it i s not very common (e.g. iLibeskind et alJ 
200^: Kang^t aljJ^O^ ZentoeretjtL J200^ : iLi & Helmil 2008^ 
Libeskind et al.ll2009l : lDeason et alj|201llk 



4 STAR FORMATION HISTORIES 

4.1 A comparison to Local Group dwarf galaxies 

The resolution for which the stellar ages can be determined depends 
on the populations that can be used as tracers. Overall, the resolu- 
tion decreases with increasing age and ranges from a few Myr at 
the youngest end (ages up to 1 Gyr), to several Gyrs for stars older 
than a few Gyrs. In this work we consider three different age bins 
that can be well separated in a CMD analysis: an old population 
(> 10 Gyr), an intermediate population (1-10 Gyr) an d a young 
population, < 1 Gyr (see also lTolstov, Hill & TosillSlO^ . 

In the top left panel of Figure we show the percentage of 
satellite galaxies that contain observable populations (defined as 
> 1 per cent of the total mass) in each of the age bins in the simu- 
lations Aq-A to F, where red asterisks correspond to the old, blue 
diamonds to the intermediate and green triangles to the young pop- 
ulations of stars. 

From this panel we see that ~9 per cent of all modelled satel- 
lites do not have an old population, whereas all Milky Way satel- 
lites and all isolated dwarf galaxies observed with sufficiently deep 
colour-magnitude diagrams do contain old stars. We have checked 
explicitly that all of these subhaloes have appeared in the simula- 



tion well before a lookbacktime of 10 Gyrs, therefore the lack of old 
populations is probably related to the implementation of the semi- 
analytic prescriptions. There are at least two possible explanations 
for this difference associated to how we model star formation and 
cooling. 

Firstly, we do not allow cooling in dark matter haloes with 
Tvir < 10"^K. If, however, we relax this assumption and also al- 
low cooling in such haloes with the same efficiency as a halo of 
Tvir = lO^K, we find that a larger percentage of all galaxies (95 per 
cent of all dwarf galaxies, and 97.5 per cent of those that are satel- 
lites presently) does form an old population. But as a consequence 
the faint end of the satellite luminosity function is also much en- 
hanced (most notably at Mv > —5). A proper implementation of 
the physical processes playing a role in the formation of the first 
stars would require many extra assumptions, for instance on the 
IMF and the interplay between H 2 cooling and H 2 dissociation in 
the host haloes which we prefer not to include in our models at this 
stage. 

Another limitation of our semi-analytical model is that we do 
not represent the stochasticity of star forming regions. In our mod- 
elling, for stars to form in a disc, the total (global) surface density of 
the disc has to be above the star forming density threshold. How- 
ever, in reality star formation processes can be much more local, 
i.e. one molecular cloud can have the required density while its 
surroundings might not. 

In the left middle panel we plot the percentage of galaxies 
for which a given population is dominant (i.e. more than 50 per 
cent of the stars belong to it) as a function of My. Table sum- 
marizes the available data for the Milky Way satellites (the fir st 5 
columns are taken from Table 2 of iTolstov, Hill & Tosill200^ . A 
direct comparison can be made for galaxies in the luminosity bin 
— 12 <Mv < —8.5. Most modeled satellites in this bin are dom- 
inated by an old population, while a significant minority (^20 per 
cent) is dominated by an intermediate population of stars. This is 
very comparable to the results for the satellite galaxies around the 
Milky Way, as can be directly seen from their percentages over- 
plotted in the left middle panel of Figurej^as filled circles, where 
out of the seven galaxies in this luminosity bin only two are clearly 
dominated by an intermediate population. 

The left bottom panel of Figure shows the relative contribu- 
tions of old, intermediate and young populations averaged for all 
galaxies in a particular luminosity bin. 

While all left panels show the satellite galaxies of the main 
Aquarius galaxy, the right panels of Figure show the contribu- 
tions of the various populations for dwarf galaxies outside the main 
halo, i.e., from 400 kpc to 2 Mpc. We focus here on slightly brighter 
dwarf galaxies, with — 19 < My < —8.5, as these may be ob- 
servable with current instrumentation. Most of these systems are 
the main galaxy within their dark matter halo, which are the only 
galaxies fed by cooling. The isolated dwarfs are mostly dominated 
by intermediate age instead of old populations. Also, a young pop- 
ulation of stars is present much more frequently (compare the two 
top panels). Qualitatively, it is clear that our model produces an 
age-density relation in agreement with the observations of Local 
Group satellites. Gas-poor and older galaxies are found in over- 
dense regions near bigger galaxies, while the star forming galaxies 
are found in isolation. This result is more clearly illustrated in Fig- 
ure]^ where we show the location of —19 < Mv < —8.5 galaxies 
overlaid on the dark matter distribution for the Aquarius A simula- 
tion. (We have checked that the results shown in Figure are not 
influenced by projection effects.) The galaxies dominated by an old 
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Figure 9. The percentage of galaxies in each luminosity bin which have 
> 1 per cent (top panels) old (red asterisks), intermediate (blue diamonds) 
and/or young (green triangles) populations. In the middle panels the lines 
with similar symbols indicate whether the old, intermediate or young pop- 
ulations are the dominant one in the galaxy (that is, more than 50 per cent 
of their stars originate from this epoch). In the bottom panels we show the 
fractions of the old, intermediate and young stars averaged over all galax- 
ies in the luminosity bin. In all left panels the satellite galaxies, within 280 
kpc of the main halo are displayed, in all right panels the isolated galaxies. 
Overplotted in the left middle panel are the percentages of these Milky Way 
satellites with — 12 <Mv < — 8.5 dominated by either old or intermediate 
population as filled red and blue circles respectively for a direct comparison. 



stellar population show a very different distribution as function of 
the density of their environment than the intermediate age galaxies. 

Some galaxies outside the main haloes are dominated by old 
populations. In all of these some stars are formed at intermediate 
ages as well, but there are fewer of these than old stars. All show 
a bursty star formation at intermediate age. Some have never been 
satellites in their life, so this behaviour is entirely due to their own 
internal feedback processes. In the Local Group we also find a few 
examples of old and passively evolving small galaxies, like the Ce- 
tus and Tucana dwarf spheroidals, with no clear association to the 
Milky Way or M31, although a past association can not be ruled 
out. The existence, or non-existence, of truly isolated quenched 
dwarf galaxies in observational data will provide important con- 
straints on the modelling of the star formation threshold for dwarf 
galaxies. 

A more quantitative comparison between isolated dwarfs in 
our models and in the nearby Universe is difficult. Many of the 




Figure 10. Mass-weighted age as a function of their distance of the satel- 
lite to the main galaxy. The black circles show all Aquarius satellites with 
Mv < —5 and the grey line represents the mean in this sample. Red 
asterisks represent mass-weighted ages for the Milky Way satellites from 
lOrban et al.l j2008ll . Note that these are derived from HST observations 
and, because of the limited field of view and population gradients known 
in dwarf galaxies, might he biased towards lower mass-weighted ages. The 
two blue triangles are the Magellanic Clouds. The smaller red asterisk rep- 
resenting a lower limit is the Sagittarius galaxy, which has been severely 
stripped and therefore will have lost preferentially part of its older popula- 
tion. 



well studied galaxies outside the Milky Way are still associated 
with its nearest neighbour, Andromeda, or with the Local Group 
environment as a whole. Beyond the Local Group, at ~ 1.3 Mpc, 
the observational CMDs are much harder to interpret, since only the 
giant branch is bright enough to be resolved. We have made a ten- 
tative comparison with the ACS Nearby Galaxy Survey (ANGST) 
survey, a project measuring the star formation histories of galax- 
ies outside the Local Group out to 4 Mp c in a systematic manner 
dPalcanton et alj2009l : IWeisz et al.l20T 3). Only three galaxies have 
reliable measurements of the dominant population in their sample 
from 1.3 to 2 Mpc. Of these three galaxies, two are dominated 
by an intermediate population and Just one by an old population 
dWeisz et alj|201 1[) . 

Figure [T^ shows that there is a clear trend in mass-weighted 
age with distance from the host. Systems dominated by interme- 
diate populations are found preferentially in the outskirts of the 
host halo while older systems are found closer in. We see a sim- 
ilar trend in the Milky Way satellites, where the galaxies domi- 
nated by intermediate populations are located at greater distances 
than ~ 100 kpc as indicated in Table (with the exception of the 
Magellanic Clouds and the heavily disrupted Sagittarius galaxy). 
Note that the mass-weighted ages for the Milky Way satellites 
plotted here are derived from HST observations, with a relative 
small field of view. These mass-weighted ages will therefore be 
biased towards the more concentrated younger population if a ra- 
dial age gradient is present in the galaxy and metallicity and age 
gradi ents have been observed in many dwarf spheroidal systems 
(e.g.. lHarbeck et alj200ll : lTblstov et alj2004[B^attaglia et alj2006l : 
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Name 


Mv 


Detected Pops 


Dom. Pop 


ref 






old 


im 


yng 






LMC 


-18.1 


yes 


yes 


yes 


old 


[1] 


SMC 


-16.2 


yes 


yes 


yes 


im 


[2] 


Sagittarius 


-13.4 


yes 


yes 


no 


7 


[3][4] 


Fornax 


-13.2 


yes 


yes 


no 


im 


[5][6] 


Leo I 


-11.9 


yes 


yes 


no 


im 


[3][7][8] 


Sculptor 


-11.1 


yes 


yes 


no 


old 


[3][9] 


Leo II 


-9.6 


yes 


yes 


no 


7 


[3][7][10] 


Sextans 


-9.5 


yes 


yes 


no 


old 


[11][12][13] 


Carina 


-9.3 


yes 


yes 


no 


im 


[3][7][14] 


Ursa Minor 


-8.9 


yes 


yes 


no 


old 


[3][7][15] 


Draco 


oo 

CO 


yes 


yes 


no 


old 


[3][16] 


CVnI 


-8.6 


yes 


yes 


no 


old 


[17] 


Hercules 


-6.6 


7 


7 


7 


7 




Bootes 


-6.3 


yes 


yes 


no 


old 


[17] 


UmAI 


-5.5 


7 


7 


7 


7 




Leo IV 


-5.0 


7 


7 


7 


7 




CVn II 


-4.9 


7 


7 


7 


7 




Leo V 


-4.3 


7 


7 


7 


7 




UmAII 


-4.2 


yes 


yes 


no 


old 


[17] 


Com Ber 


-4.1 


7 


7 


7 


7 




Boo II 


-2.7 


7 


7 


7 


7 




Willman 1 


-2.7 


7 


7 


7 


7 




Segue 1 


-1.5 


7 


7 


7 


7 





Table 2. Populations in Milky Way satellite galaxies. References are: [1] 
iHarris & Zaiitskvl j2QQ9h (but note that the intemiediate and old stars are 



almost contributing equallv in the LMC SFH), [21 iHarris & Zaritskv 


2004 


). [31 iDolohinl J2OO2I). [41 IBellazzini. Ferraro & Buonannc 


1999 


). 151 IColeman & de Jond J2008h. 161 


Gallart et al.l J2005l). 171 


Hernandez. Gilmore & Valls-Gabaud J2OOOI). 


81 IGallart et al.l h999li. 


[9] Jde Boer et al.ll2012l). [101 iGullieuszik et al 


20081). [Ill Mateo etal. 


h992h. [121 iBellazzini. Ferraro & Pancinol 


200 


1 ). 1131 iLee et al.l 


2003). 


[14] 

[16] 


Hurlev-Keller. Mateo & Nemed Jl99f 


). [151 ICarrera et alj ' 


2002). 


Anaricio. Carrera & Martmez-Deleado 


teOOll). [171 IdeJonsetal. 


|2008 


). 







[Bernard et al .120081:1 |d£^oere^ J2012h as well as along the Sagit- 
tarius stream llBellazzini et al.ll2006i e.g.l. 



4.2 The physics shaping the modeled star formation histories 

Figure fTTI illustrates the variety of star formation histories found 
for the different satellites in Aquarius halo B. The 18 most lumi- 
nous satellites are shown down to a magnitude of Mv=-7.9. The 
star formation rates have heen normalized in this figure to the high- 
est peaks, hut the absolute values range approximately from 0.2 
M 0 /yr in the top row to 0.005 M 0 /yr in the bottom row. As can 
be seen from the figure, most stars in our satellites are made af- 
ter reionization (at 2 = 11.5 in our models), but before the time 
that the galaxy fell into the main halo (and thus became a satel- 
lite) as indicated by the vertical dashed blue lines. Taking all the 
stars in satellites brighter than Mv = — 5 and within 280 kpc in 
Aquarius haloes A-F, we find that 99.7 per cent of the stars formed 
after reionization, and 96.8 per cent before the satellite fell into the 
halo. The small satellites make a larger percentage of their stars be- 
fore reionization. This is because they have just a few small bursts 
of star formation each of which contributes a significant fraction 
of the stars. In the larger galaxies the star formation rate is much 
higher and it extends over longer periods so that the initial episode 



of star formation (i.e. before reionization) is not important in a rel- 
ative sense. On the other hand the quenching of star formation after 
infall ensures that a relatively small amount of stars will be made 
afterwards. 

Infall onto a host system in our model quenches star forma- 
tion, because we assume the hot halo of a galaxy is stripped as soon 
as the galaxy becomes a satellite and the hot component is added 
to that of the host. Some m odels follow ram pressure stripping pro- 
cesses mo r e gradu a lly (e.g. Font_et_a]j200^ Weinmann et al.l2010l: 

I Guo et al] 1201 Ibl: iN ichols & Bland-Flawthoml 1201 ih . iFont et alJ 
d2008h and Guo et alT ~i 201 Ibll show that their implementation 
does not make a significant difference for small systems as con- 
sidered here. Some fully hydrodynamical simulations have also 
shown that the remaining gas in the satellite galaxies is (almost) 
completely stripped onc e a galaxy falls into the main halo (e.g. 
lOkamoto & Frenl3l2009l : lokamoto et al.ll2010h and that the cool- 
ing of gas onto sa tellites seems to b e small and important only for 
massive systems dSaro et an i2010h . In our model, a satellite can 
still make some stars after infall, but only as far as its cold gas 
reservoir allows it. We find that the amount of cold gas still present 
in the modeled satellite galaxies, is in almost all cases larger than 
is observed within satelli te galaxies in the Milky Way halo (e.g. 
iGrcevich & Putm^l2009h . Flowever, because this gas is below the 
star formation threshold, it does not form any stars. It is quite pos- 
sible that ram-pressure stripping of cold gas might have removed 
some of this remaining gas from the observed satellite galaxies. 
Additionally, a more realistic model of star formation, not based 
on a global density threshold, but following the local density of the 
gas, might also alter the amount of gas left in these galaxies. The 
fact that all the classical dwarf spheroidals have a detected inter- 
mediate population, whereas this population is not always present 
in the models (see top left panel of Figure |9} hints that our model 
might shut off star formation too soon after infall. 

A change in the value of the threshold density for star for- 
mation, equivalent to a more stochastic implementation to reflect 
molecular cloud physics, will not lead to continuous star forma- 
tion histories for most satellites, since the density of cold gas is 
too far below the threshold most of the time. Flowever, there are 
some cases where the star formation threshold is barely not met, 
and a small change in the value of the threshold, could significantly 
change the star formation history. The most bursty star formation 
histories are found in the galaxies which are constantly close, but 
mostly below, the star formation threshold. 



5 MILKY WAY DWARF ANALOGS 
5.1 Star formation histories 

In Figure [T^ we show some model satellite galaxies that are com- 
parable in either brightness, average metallicity and star formation 
histories, or several of these properties, to the Carina, Sculptor and 
Fornax dwarf spheroidals. For each galaxy we show three examples 
to illustrate the scatter in these properties. 

In the bottom panels of Figure d we plot the observed 
star formation history for each dwarf spheroidal. Even though 
it is possible to observe all three galaxies down to their main 
sequence turnoff magnitudes, large differences are found in de- 
rived star formation rates due to incomplete spatial coverage or 
the s et of model isochrones used (see for example iGallart et alj 
l2005h . We have ch osen to show here t h e star formation histo- 
ries as obtained by IColeman & de Jon3 J2008h for Fornax and 
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Figure 11. The relative SFRs as a function of lookback time for the most luminous satellite galaxies in Aquarius B (black solid lines). We have indicated the 
infall time (blue vertical dotted line). In the top left comer of each panel, the absolute V-band magnitude of that particular galaxy is given and the scaling factor 
by which the star fomiation rates have been multiplied to obtain the normalized values. 



by Ide Boer et alj j2012l) for Sculptor, because these are the most 
recent analyses, derived from photometry which covers (nearly) 
the entire galaxy. For Fornax, we also s how, w it h a d ashed 
line, the star formation as compiled by iGrebell dl998h . The 
star formation history of Carina ma tches the one derived by 
[Hurlev-Keller, Mateo & Nemed dl998h . with its characteristic three 
peaks. Other studies claim slightly different star formation his - 
tories llDolDhinll2002l : [Hernandez, (Tilmore & Valls-GabaudllrOOOh . 
Of course, one needs to bear in mind that the error bars on age are 
quite large for the observations (generally in the order of Gyrs for 
the regime we are most interested in), so the studies are not neces- 
sarily inconsistent. 

For Carina, we show several candidates selected on luminosity 
(—10.5 <Mv < —8.5) and their distinct bursty star formation with 
a majority of the stars formed at intermediate age. Amongst all our 
Aquarius simulations, we only find six candidates that have a dom- 
inant intermediate population in the corresponding luminosity bin. 
In Figure [T 2 I we show three candidates which have a very bursty 
star formation history. It is interesting that also bursty Carina-like 
galaxies are produced. In our models Carina’s very bursty star for- 
mation history occurs before it becomes a satellite and is the result 
of an interplay between gas density and star formation threshold, 
and is clearl y not due to tidal in teractions with the host (such as 
suggested bv IPasetto et alj|201lh . All of these Carina-analogs fall 
into the main halo quite late. As discussed before in Section |4^ 
most of these model galaxies find themselves close to the threshold 
of cold gas needed to make stars, and are only above the threshold 
occasionally. This explains the bursty nature of these star formation 
histories in our model. We have explicitly checked their merging 



history, but found no significant merging events causing or preced- 
ing any of the bursts. The Carina model shown in the bottom panel 
of Figure [T2I depicts a star formation history that matches best the 
observed general properties of Carina, i.e. also matches the metal- 
licity. Also this galaxy has a strong intermediate age burst of star 
formation, but a lower peak at a young age. 

For Sculptor we found many analogs for the star formation 
history, since a great majority of our model satellite galaxies are 
dominated by an old population. However, most candidates of com- 
parable luminosity have a slighter higher metallicity than Sculptor. 
Three typical examples of comparable luminosity, and also average 
metallicity (the largest offset is 0.5 dex), are shown in the middle 
column of Figurell2l 

The Fornax candidate model shown on the bottom panel of 
the right column of Fi gurell2lmatches very well t he observed star 
formation history froni lColeman & de long! (I2OO8I) and the general 
properties of the galaxy. On the other hand the models shown in the 
top and middle panels have star formation histories which match 
better older observations of Fornax in which the pea k of star forma - 
tion occurs at an older age, such as that compiled by iGrebeil d 1 99 8h . 
All examples shown here were chosen to have a comparable lumi- 
nosity and average metallicity (the largest offset in metallicity is 
0.2 dex). 

5.2 Metallicity distributions 

Figure shows the metallicity distributions for the same set 
of model galaxies chosen to be relatively close analogs of the 
Milky Way dwarf spheroidals Carina, Sculptor and Fornax shown 
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Figure 12. Star formation histories for models which resemble the three 
classical dwaif galaxies Carina, Sculptor and Fornax (top panels). The star 
formation histories are shown both in the original binning of the simulation 
(grey lines) and rebinned to bins of 1 Gyr, to provide an easier comparison 
to the observations. In both cases the star formation histories are normal- 
ized to their peak value. Overwritten within the panels are the metallicity, 
luminosity and the scaling factor by which the model star formation rates 
have been multiplied to get the normalized values shown here. Division by 
this scaling factor thus returns the model star formation rate values, for the 
non-rebinned model (grey lines) in M 0 /yr. We have also indicated the infall 
time of the satellite as a blue vertical dotted line. In the bottom panels we 
show observational star formation histories for the same galaxies (see text 
for references). 
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Figure 13. Metallicity distributions for the models displayed in Figure 1121 
which resemble the three classical dwarf galaxies: Carina, Sculptor and For- 
nax (top three panels). In the bottom panels we show observational metal- 
licity distributions for the same galaxies (see text for references) in [Fe/H] 
from CaT samples (black solid lines), or con'ected for a global relation to 
[Mg/H] values (grey solid line). 



discussed in Section|4] Since the first peak of star formation is very 
small, it will only enrich the galaxy hy a small amount leading to 
the formation of very low-metallicity stars in the second hurst. 

Overall, Figure [TSishows that the form of the distribution can 
vary significantly, depending on the exact star formation history of 
the system. 



in Figure [T2] The bottom panels show the observed metallic- 
ity distributions for these dwarf spheroidals as taken from the 
DART data sets. The black histograms correspond to the observed 
[Fe/H] distribution as de rived from the Ca II triplet calibration of 
IStarkenburg et al.l J2OIOI) . The grey histograms give the metallicity 
distributions corrected for [Mg/H] using the simplified relation of 
[Mg/Fe] discussed in Section [331 We regard this as a more direct 
comparison to our model, which assumed instantaneous recycling. 

In general, the metallicity distribution functions of the model 
dwarfs are narrower, although a typical error on the observed metal- 
licities is 0.2 dex, which will broaden the distribution. There is also 
a lack of extremely metal-poor stars in most metallicity distribu- 
tions of our models. In the case of the Carina model shown in the 
top left panel however, too many very metal-poor stars are pro- 
duced and the global metallicity is too low. This also is probably 
driven by our crude modelling of the first generations of stars, as 



6 DISCUSSION 

6.1 Comparison to the galform code 



The Durham semi-analytic model of gal axy formation, GALFORM 
(e.g. lCole et alj|200(il : [Bower et al.ll2006l). has al so been run on the 
same s ix Aquarius haloes by Cooper et al.l (I2OIOI) and bv iFont et alj 
( I2OI ll) . The latter implemented a novel and detailed treatment 
of both local and global reionization and compared the luminos- 
ity functions of their mode l and the satellite ejection model of 
lUi. De Lucia & Helmil ( |201(]|) . This is identical to the one presented 
here in the bottom panel of Figure |3 since in the comparison no or- 
phan galaxies are considered and ou r stellar stripping mechanism 
does not significantly alter its shape. [Font et aP (1201 1[) concluded 
that both models are quite similar in methodology as well as results, 
including our treatment of reionization which is much less sophisti- 
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cated as theirs. However, we note that a detailed treatment of local 
reionization could possibly affect also other ob servables such as 
the rad ial distribution of satellites, as shown by lOcvirk & Aub^ 

( l201lh . 

Our model and those resulting from GALFORM show some in- 
teresting differences however. For example the total stellar mass 
in the main galaxies are generally lower in the GALFORM models, 
ranging from 7 x 10® Mq to 1 x 10^^ Mq (compare Table 1 from 
iFont et aP 1201 ill and Figure [^in this work). The number of satel- 
lites around the main halo are very similar, although there clearly 
are fewer very luminous satellites (brighter than My = —15) 
present in the GALFORM models. This seems to be due to their 
slightly stronger feedback in this regime. In general the GALFORM 
satellite luminosity fu nction is slightly steeper, but consistent with 
IKoposov et alJ l2008l ) and our results, in particular for luminosi- 
ties fainter than My = —15. Their luminosity-metallicity relation 
flattens out significantly at the fainter end (see their Figure 5, right 
panel), which is an interesting result directly related to their satu- 
rated feedback scheme, which is clearly distinct from the feedback 
scheme used in this work. Unfortunately, the errors on the available 
data are too large to test this prediction. 



6.2 A stellar mass-halo mass relation: comparison to 

abnndance matching and hydrodynamical simulations 

An interesting issue in the context of the ACDM paradigm is 
which galaxies reside in which haloes. Several groups have quan- 
tified the relationship between stellar mass and halo mass us- 
ing N-body simulations combined with ( s emi-analytical) mod- 



els o r observations (e .g. Frenk et al.| 19881: Navarro & Stei nmetj 
200d: Dekel & 2003 : Yang. Mo & va n den Bosch 20031: 



Wang et a lj|2006l: Shankar e talJl2006L Yang. M o & van den Bosch 



Conroy & Wechslei 20091 :" Moster et alT 2010l; Guo et al 



2008 



20101) . This comparison is made under the assumption that the dark 



matter halo mass and the stellar mass of a central galaxy follow 
a monotonic relation, taking i nto account some amount of scatter. 
For example IGuo et af] 1 I 2 OIOI) combined the abundance of galax- 
ies with stell M masses in the r ange 10® M© < M* < 10^® Mq 
from SDSS jUi & Whitel |2009|) . and the dark matter halo abun- 
dances from the Millenniu m Simulations dSoringel et al.l l2005l : 
iBovlan-Kolchin et alj|200d ). The relation derived from this anal- 
ysis predicts a dark matter mass for a Milky Way-type galaxy of 
~ 2 X 10^®Mq, which is a factor of two larger than the favoured 
result from our semi-analytic model, but within the current obser- 
vati onal constrain ts. 

ISawala et all ) l201lh noticed that current hydrodynamical sim- 
ulations of isolated dwarf galaxies produce an order of magnitude 
larger stellar masses than the expected relation from IGuo et alj 
l l2010h extrapolated t o the lower halo mass end. In Figure IT? we 
show the relation from iGuo et al.ld2010l) an d iMoster et al.l (1201 01) as 
well as the results predicted by our model for all the central galax- 
ies in the Aquarius simulation^] and the hydrody n amical simu- 
lations from| Pelup essv. van de r Werf & I ckj j 2004]): IStinson et alj 
| |2()07 ): Govemato et 1 201(]lk1sawala et al.l l 201ll) . as compiled 
bv ISawala et al, < 2011 . see also their Figure 4.). 



We see a strong correlation between total stellar mass and 
dark matter halo mass (taken as their virial mass at z = 0, as in 
IGuo et ani2Ql(]|) for the central galaxies in our model, as shown 



^ We have used a distance limit of 2.5 Mpc from the centre of the box to 
stay well within the high-resolution region of our simulations. 




1q9 ^qIO ^q11 ^q12 ^q13 

DM halo mass [M©] 



Figure 14. The relation between stellar masses and dark matter halo 
masses for all central haloes within 2.5 Mpc from the main halo in 
the Aquarius simulations (black filled circles). The red symbols are 
results from hydrodynamical simulations, as compiled in [Sawala et alj 
i201lh. fromlPelupessy. van der Wer f & lcl3 l l2QQ4l) : IStinson et alJj2QQ7h : 
iGovemato et affboiOh lSawal^et al.l j2Ql ll) . The blue asterisks are galax- 
ies from an SDSS sample for which stellar masses were derived from spec- 
tr oscopy and dark mat t er viri al halo masses from weak gravitational leasing 
bv iMandelbaum et al.U2QQ^ : the error bars give the 95 per cent confidence 
intervals. The black solid line re presents the stellar mass-halo mass relation 
as derived bv iGuo et al.l ll2Qld) . constrained by SDSS DR7 data. The ex- 
trapolation into the low-mass regime is indicated by the black dashed line, 
whereas the light grey area shows the maximum dispe rsion, <T]og-M.> = 0. 2 
dex. The stellar mass-halo mass relation derived by iMoster et alj j2Q10h . 
including a scatter of <TiogM* = 0-f^ (dark grey area), is shown as a 
black dashed-dotted line and its extrapolation into the low-mass regime is 
given by the dotted black line. 



in Figure Q4] Howe ver, the relation p redic ted from our model and 
the one predicted bv iGuo et al.l bOlOh an d iMoster et al.l d2010l) are 
clearly offset: at a given luminosity our model galaxies reside in 
smaller dark matter haloes. T he difference is gre atest for the lower 
mass ohjects, for which the IGuo et aO (120101) and iMoster et alj 
d2010l) relations are extrapolations (i.e. these scales are poorly con- 
strained by the data used by these authors; and are indicated by 
the dashed and dotted lines in Figure [T4l l. but also at the regime 
of the main central galaxies (our ‘Milky Way’s) there is a signifi- 
cant difference (with the exception of the main galaxy in Aquarius 
E). Ge nerally, our results are in agreement with Mandelbaum et al] 
d2006h . who used SDSS to derive stellar masses from spectroscopy 
and virial halo masses from weak gravitat ional leasing. Inte rest- 
ingly, while bein g offset from the results of iGuo et al.l d2011bl) and 
the relation from iGuo et al. I J 20 IT at the lower masses, our model 
predictions are also consistent with most hydrodynamical simula- 
tions showm 

Both IMoster et alJ d201(]|) and IGuo et al] d2O10l) demonstrate 
that the relation for the most likely stellar mass within a given dark 
matter halo changes with different assumptions about scatter. The 
scatter wit hin the Aquarius halo es is larger than th e predictions 
from either lOuQ et al.ld2010l) and IMoster et al.l d2010l) . certainly at 
the lowest masses, raising the question of whether environmental 
stochastic effects (such as those associated to the time of formation, 
the mass at reionization, the cold gas density and the star formation 
threshold) may play a more important ro l e for lo w mas s systems. 

Re cent work of IVera-Ciro et al.l d2013l) and IWang et ^ 
d2012bl) show that the assumption of a lower mass for the Milky 
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Way aro und ~ 8 x 10^^ will alleviate significantly th e problem 
raised bv iBovlan-Kolchin. Bullock & KaplinghaJ 1201 il l that satel- 
lite systems around Milky Way-sized haloes would be too dense 
to host the known Milky Wa y dwarf galaxy satelli tes. Using the 
model presented in this paper IWra-Ciro et al.ld2013l) subsequently 
show that good agreement can be found for the relation between 
luminosity to total mass of the modeled satellites in comparison 
with observational constraints (in disagreement with the Vmax-Mv 
dichotomy between modeled and observed satellites as shown in 
Figure 6. of lBoylan-Kolchin, Bullock & Kaplinghai 20121 , based on 

abundance matching). 

Using an ext ension of the lPe Lucia & Blaizod j2007t ) model, 
IGuo etakfeoilbl) obt ain a g ood correspondence with the relation 
derived by IGuo et all 1201 OP for SDSS galaxies and consistency 
with galaxy luminosity functions over a large range in magnitu des 
in various bands and with the extension of the lGuo et al.l 1 I 2 OIOI) re- 
lation for lower mass galaxies. They found that in order to reach this 
good correspondance they needed to adopt a supernova feedback 
ejection efficiency which depends on the circular velocity of the 
underlying dark matter halo with an exponent = —3.5, whereas 
in our model /3i = —2. A careful comparison of both models on 
larger cosmological scales is beyond the scope of this paper and 
also not feasible with the set of simulations used, but we plan to 
address this issue in future work. 



7 CONCLUSIONS 



In this work we have confirmed that a semi-analytic model of 
galaxy formation applied to a high resolution cosmological N-body 
simulation is able to match observed relations on the scale of the 
Milky Way and its satellites simultaneously. We have compared the 
results for the luminosity function, luminosity-metallicity relation 
and radial and spatial distribution to observations of dwarf galaxies 
around the Milky Way and investigated the age-density relation and 
star formation histories within the model. 

We find that the same Aquarius halo in which we found the 
closest Milky Way-like galaxy analogue, halo B, also provides us 
with the best matching satellite system in terms of its luminosity 
function. This galaxy does not resemble the Milky Way in every 
respect. For example the radial distribution of satellites within this 
particular simulation is more centrally concentrated than the ob- 
served distribution of satellites around the Milky Way and it does 
not show the same degree of spatial flattening. It also does not host 
a galaxy with a luminosity comparable to the LMC or SMC. Of 
course, the Milky Way, although in many ways the best studied 
galaxy we know, is j ust one example and it might be u nusual in 
several respects fe.g.. lFlvnn et al.ll200 ^ IGuo et al.ll201 lah . 

We find a clear relation in our models between the num- 
ber of bright satellites and the host dark halo mass. With our 
current feedback and reionization prescriptions, our best Milky 
Way analogue has a dark halo close to 8 xlO^^ Mp, in agree- 
ment with the lower estimates from ob servations feattaglia et all 
l2005l : [Smith et al II 20071 : IXue^talj [200^. but a factor of two lower 
than the best estimates from lUi & Whitel 120081) : IGuo et al.l (1201 Ot) . 
Additionally, over the mass range probed by the Aquarius sim- 
ulations, we find a different relation between dark matter halo 
m ass and ste l lar m ass o f the central galaxi es than that derived 
by Guo et all ll2010t) and iMoster et al.l 120101) . and reproduced by 
IG^uo et ^ N2011b|) . 

Based on their star formation histories, we find model 
satellites analogous to the Sculptor, Fornax and Carina dwarf 



spheroidals, although none of the model galaxies provides a match 
of all observable properties. However, the metallicity distributions 
for these galaxies are generally too narrow compared to the obser- 
vations and they lack an (extremely) metal-poor population. It is 
unclear at the moment whether this can be completely ascribed to 
the lack of a detailed prescription of the chemical evolution of dif- 
ferent elements and to the adoption of the instantaneous recycling 
approximation. This topic will be the subject of further work. Also 
our model does not allow cooling via H 2 in haloes below the Hy- 
drogen atomic cooling limit, and hence does not provide a fully 
physical model for the formation of the first stars. Another short- 
coming of our current model is that our modeled satellite galaxies 
have too much cold gas compared to observations. This could per- 
haps be (partly) solved by including ram-pressure stripping of the 
cold gas when an object becomes a satellite. 

However, various predictions can be made from our model that 
are expected to be independent of these shortcomings. We predict 
the ratio between galaxies dominated by old- or intermediate popu- 
lations of stars to be close to 1:2 beyond the Local Group. We also 
expect a large majority of the satellite galaxies that are dominated 
by intermediate age stellar populations, to have fallen into the main 
halo relatively late, including galaxies with a bursty star formation 
like Carina. Generally, a very small percentage of stars in a satel- 
lite is formed after its infall-time, due to stripping of the satellite’s 
hot halo which prevents further gas from cooling and forming stars. 
However, some galaxies within the model have stopped star forma- 
tion several Gyrs before their infall, due to internal processes, such 
as supernova feedback, which have expelled their gas. Our model 
thus predicts the last major star formation event within a (classical) 
satellite to rather provide a lower limit on the time elapsed since 
infall. 

The brighter galaxies amongst the ultra-faint satellites (e.g 
—8.5 <Mv < —5), which could be well resolved in our model, 
have formed a larger percentage of their stars in a single burst, 
because feedback has a larger impact and prevents a continuous 
mode of star formation. These model ultra-faint galaxies are gener- 
ally older than the more luminous counterparts, and hence contain 
a higher fraction of stars formed around, and even before, the epoch 
of reionization. 
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APPENDIX A: IMPLEMENTATION OF STELLAR 
STRIPPING AND TIDAL DISRUPTION 

AI Stellar stripping 

In our simulations, the dark matter subhalo of a satellite may be 
so heavily tidally stripped that also its stellar component should 
be affected. Here we implement an approach to ensure that no 
very extended satellite galaxies reside within much smaller, heavily 
stripped dark matter subhaloes. This approach is thus only used on 
satellites which still have a dark matter component. 

The half-mass radius of an exponential disc is related to its 

scale-length through R\ = 1.67 x R^. At each time-step we now 

2 

compare the half-mass radius of the stars (and cold gas), R\ , to the 

2 

half-mass radius of the dark matter subhalo, The value of 

is measured by the SUBFIND algorithm dSpringel et al.l200lh 

2 

directly using the dark matter particles belonging to the satellite. If 
through tidal stripping in the simulations 7??^ < R\ , we remove 

the stars (and corresponding cold gas) up to R^^ . This implies 

2 

that after stripping the disc now has a mass given by 

rY^ 

Md-new = 2Md[l- (l-f-^)e i ]. (Al) 

The updated exponential disc, consisting of the leftover cold gas, 
has a half-mass radius set by Rl — R^^ and the new scale 

2 --new 5 

radius for the disc is therefore i?d-new = 7(?^/1.67. 

2 

The cold gas and stars stripped from the satellite, are trans- 
ferred respectively to the cold gas disc and spheroidal stellar com- 
ponent, which includes the bulge and stellar halo, of the main 
galaxy. 

Although the assumption that 50 per cent of the stars should 
at least be contained within the dark matter half-mass radius seems 
like a natural first guess, there are very little direct observational 
constraints on the exact ratio. However there is general agreement 
that the dark matter halo indeed has to be severely stripped before 
stars are affected. For instance, the scale length of the Milky Way 
stellar thin disc is ~ 3 kpc, whereas its dark matter half-mass ra- 
dius is thought to be in the order of 100 kpc. For a more in dept 
discussion of tidal stripping in disc galaxies, and the influence of 
parameters as galaxy mass, orbit and disc inclination from n umeri- 
cal simulations we refer the reader to lVillalobos et al.ld2012h . 



A2 Tidal disruption 

In our tree files a dark matter subhalo is ‘lost’ when it is stripped 
down to less than 20 particles (corresponding to a dark matter mass 
of ~ 2.7 X 1 0^ Mpi in the Aquarius A halo simulated at level 2 
ISpringel et al.ll2008ah . We refer to these galaxies as ‘orphan’. The 
question arises: when their dark matter subhalo disappears what 
should be their fate. 

We assume that the satellite galaxies residing at the cen- 
tres of the fully (below the resolution limit) stripped subhaloes 
will merge with the c entral galaxy after a dyna mical friction time 
scale, as explained by |Pe Lucia & Helmil l l2008l) . However, for low 
mass galaxies this time scale will be longer than the Hubble time. 
In some pr evious work on semi-analytical models of Milky Way 
haloes (e.g. lUi et alJl200lj : lLi. De Lucia & Helmill20ld : lForit et alj 
l201lh these orphan satellites galaxies were not taken into account, 
since they were believed to get completely stripped. However this 



needs to be checked explicitly and in some cases we can expect the 
galaxies to survive the tidal forces. 

Here we improve further the code, by comparing the average 
mass density of orphan satellites to that of their host system and 
determine whether they will survive. We do this only for those or- 
phan galaxies that would otherwise survive, so will not merge on a 
dynamical friction time-scale with the central galaxy before z = 0. 
This means we treat the more and less massive orphan galaxies 
via different mechanisms. The more massive orphans are dragged 
in and merge with the central at the dynamical friction time-scale, 
the less massive might disrupt within the halo (or not, depending on 
their density). Differences between the two implementations appear 
firstly in the time-scale of merging (which is set at the dynamical 
friction time-scale for the massive orphans and just taken as the 
first snapshot at which they become orphans for the less massive 
galaxies) and secondly their impact on the host galaxy: only within 
the merging prescription black hole growth and the impact of major 
mergers are modeled (which can induce starbursts and morpholog- 
ical changes). Tests showed that if all orphans be treated such to 
disrupt within the halo, the final properties of the host galaxies can 
be affected. In most cases, the changes are small and mainly affect 
the bulge-to-disc ratio and the mass of the central black hole. 

For the disruption mechanism we compare the density of these 
satellites, at the time they become a orphan galaxy, to the den- 
sity of the host galaxy at pericenter. In the case that the den- 
sity of the satellite is not sufficiently high compared to the en- 
vironment, they are (assumed to be) completely disrupted. A 
similar tidal disruption mecha nism was already followed in sev- 
eral semi-analytical codes (e.g. lMo naco. Fontanot & Taffonill^OTl : 



Henrigues. Bertone & ThomasI l2008l : Guo et alj|2011bll . The ap- 
proach describe d here is in essen ce most comparable to the im- 
plementation of IGuo et alj ( l2011hh . but there are some significant 
differences, mo st importantly we use a Navarr o, Frenk and White 
profile (NFW; Navarro , Frenk & Whitel 1 1996b for the host halo 
whereas IGuo et al.l d201 Ibl) assume an isothermal sphere. 

To establish whether a oiphan satellites will survive or not, 
we need to estimate its density at the pericenter of its orbit. This 
distance may be derived numerically from: 

1 

Etot ~ h (E>NFw(77peri), (A2) 

2 Rperi 

here the total energy, Etot, and angular momentum, L, are those 
from the last recorded time. 

Subsequently the density of the satellite, (psat) is compared to 
the average density of the host, (phoat), at the pericenter distance: 

l^halo (^Rperi ) -f A7* -t- Adgas 



{Phost) — 



R^ 

pern 



(A3) 



Since the satellite’s dark matter content has fallen below the reso- 
lution limit of the N-body simulation it is not possible to evaluate 
its density precisely. Here we adopt two approaches to bracket the 
true density of the satellite. 

• To obtain an upper limit to the satellite’s density ({psat)), we 

measure the average dark matter density within in the last 

2 

snapshot the satellite was still detected within the simulation. To 
compute the satellite’s average density, we also take into account a 
fraction of the stellar and cold gas mass of the satellite galaxy. This 
fraction of baryonic mass within 77?'’^ is determined as: 



/ = 



Md(77?M) 

-mT- 



77?’ 









(A4) 
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(Psat) — /^DM -13 ■ 

• To obtain a lower limit to the satellite’s density, we assume all 
dark matter is stripped and we take as the average satellite density 
that given by the full stellar and cold gas mass of the galaxy within 
five disc scale lengths (within which 95 per cent of the stars should 
be contained). In this case: 



(Psat ) 



M* 



(5i?d,a 



0 ^' 



(A 6 ) 



We consider the satellite galaxy to be disrupted when (psat) < 
(Phost), where (paat) is given by either Eq. (ImJ or (IA61 . The cold 
gas component of the disrupted satellites is added to that of the 
main galaxy. We tested that the addition of this gas to either the 
cold or the hot gas component of the main galaxy does not sig- 
nificantly affect the properties of the main galaxies, since it might 
not be completely physical to assume that all the cold gas from the 
satellite flows into the cold gas reservoir of the host. The stars from 
the disrupted satellites are added to the spheroid component (which 
includes the stellar halo of the main galaxy). 



A3 The effects of the stripping and disruption prescriptions 

In Figure lAII we show the mass functions of present-day masses 
for the systems of subhaloes, and the mass- and luminosity func- 
tions of luminous satellites in Aq-A-2 using different modeling 
prescriptions. From the middle panel of the figure, where both the 
total subhalo and luminous subhalo mass functions are overplot- 
ted for Aq-A-2, it is clear that all massive halos at present day 
(Mbound >~ 1 O®M 0 ) will have formed stars, but that at the lower 
mass end subhalos with similar present-day masses (which can be 
affected also by stripping processes) might or might not have a lu- 
minous component. 

Figure lAll also shows the effect of the stellar stripping and 
tidal disruption prescriptions, introduced in Section lAll and IA2I 
on both the mass function and the luminosity function in Aq-A-2. 
Overall, the stripping mechanism implemented in the model has a 
small effect on the total luminosity function, only affecting a small 
percentage of the bigger satellite galaxies (Mv ~ —10). This is 
still true if we relax our assumption on the percentage of stars to be 
contained within the dark matter half-mass radius, and change this 
from 50 per cent to 20 or 80 per cent. In Aq-A-2, just 3 present-day 
galaxies with Mv < — 5 are affected by the stripping mechanism. 
Nonetheless, the stellar stripping criterion is an important addition 
for the study of individual galaxies since it prevents the unphysical 
presence of too massive galaxies within very heavily stripped dark 
matter subhaloes in the model. It also strips down galaxies which 
will become orphans eventually, making the comparison of their 
densities to that of the halo more realistic. 

The implementation of tidal disruption affects the orphan 
satellites, but not all of them, as shown in Figure I aTI 



APPENDIX B: THE EEEECTS OE NUMERICAL AND 
TIME RESOLUTION 

BI Numerical resolution 

The Aquarius Project haloes have been run at different numerical 
resolution levels and these runs show remarkable good convergence 
on the properties (mass, position and kinematics) of the simulated 
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Figure Bl. Luminosity function for Aquarius A levels 2-5 for the two ways 
of estimating the average density of a satellite for tidal disruption (see Sec- 
tion lAU . In the top panel the density estimate includes the contribution of 
the dark matter, while in the bottom panel it does not. The dashed vertical 
line and grey area in the bottom panels indicate the luminosity of satellite 
galaxies at which Aq-A-2 and Aq-A-3 start to deviate and thus the numeri- 
cal resolution be considered to affect the results. 



haloes dSpringel et al.ll^08ah . This enables us to explore the ef- 
fects of numerical resolution on our models very directly. Figure 
IB 1 [ shows the effect of numerical resolution for both tidal disruption 
approaches on the luminosity functions. The luminosity functions 
of Aq-A-2 and Aq-A-3 start to diverge significantly for My > — 5. 
This justifies the choice used throughout the article to use this mag- 
nitude limit in our comparisons to the Milky Way satellites, which 
includes many of the (brighter) ultra-faint satellites. Lower level 
resolution merger trees are only available for halo A, we can there- 
fore not perform similar convergence tests for the other simulations. 

We also choose to use the disruption prescription without dark 
matter (i.e., based on the lower limit estimate of the density of the 
satellite) in the rest of this work. Either choice would be an approxi- 
mation, due to the intrinsic uncertainties in the modelling of orphan 
galaxies. However Figure IaTI shows that the tidal disruption with- 
out dark matter shows a convergence between Aq-A-2 and Aq-A-3 
down to fainter magnitudes. At the level of My = —5 the lower 
limit approach predicts the disruption of 20 additional galaxies (see 
Figure IaTT i. while in Aq-B the difference is only 7 systems, since 
this halo has a much smaller population of oiphan galaxies. 
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Figure Al. Leftmost panel: Mass functions for all systems of subhaloes within 280 kpc of all six Aquarius haloes. Right panels: Mass (left) and luminosity 
(right) f unctions for all luminous sate llite galaxies within 280 kpc of the main galaxy in Aq-A-2. Different models are used: the default model (‘ejection’ 
model of lLi, De Lucia & Helmi|[?01Cl hlack solid line), only the galaxies with a dark matter subhalo within the default model (grey solid line), default model 
now including stellar stripping (red solid line), and default model including both shipping and tidal disruption with or without dark matter (providing an upper 
and lower limit to the luminosity function, blue and green solid lines respectively). Additionally, the mass function of the system of all subhaloes is overplotted 
as a dashed black line in the left panel. The effect of stellar stripping alone, shown here as the difference between the black and red solid line, is not visible 
in the mass function (as one would expect) and has only a small effect on the luminosity function. Notice that the masses shown are the present-day masses 
obtained by adding all bound paificles of a subhalo. This value can be affected by stripping processes and therefore be significantly reduced from the mass at 
infall. For the orphan galaxies the mass is used from the last snapshot in which they were found. 



APPENDIX C: RESOLUTION OE TIMESTEPS 

Additionally, we test the dependence of the star formation histo- 
ries obtained in the model for their dependence on the resolution of 
timesteps taken. We use Aq-A-4 for this purpose, which has a very 
large number of snapshots and can thus be rebinned into different 
time resolutions. We find that the less luminous galaxies, which of- 
ten have the more bursty star formation histories, are the most vul- 
nerable to a change in the number of time-steps taken in the simu- 
lation. By changing the time resolution in the simulation by a factor 
of two and rebinning the final output to the lowest time resolution, 
typically ~30 per cent of the stars in a dwarf galaxy will be formed 
in a different snapshot. For the galaxies with Mv < — 8.5 only, this 
shrinks to 25 per cent. In all cases, this has hardly any effect on the 
characterization of the star formation history in old, intermediate, 
or young bins as used in this paper. Only typically ~ 1 per cent of 
stars ever will change between these bins in the tests. We have also 
checked the time resolution within the semi-analytical model itself 
and found that the resolution used is sufficient; i.e. by increasing 
the number of time steps by a factor 2, changes to the results as 
used in the paper are negligible. 
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